package com.chaintor.demo.math;

import org.apache.commons.math3.analysis.interpolation.SplineInterpolator;
import org.apache.commons.math3.analysis.polynomials.PolynomialFunction;
import org.apache.commons.math3.analysis.polynomials.PolynomialSplineFunction;



public class CubicSplineCustomKeys {

    private double[] knots;
    private PolynomialSplineFunction dspsfunc;       //三次样条函数
    private double[][] pc;
    private PolynomialSplineFunction derivativeFunc; //三次样条函数的导函数

    /**
     * 自定义关键点构造函数
     * @param x
     * @param y
     */
    public CubicSplineCustomKeys(double[] x, double[] y){
        SplineInterpolator spline = new SplineInterpolator();
        this.knots = x;
        dspsfunc = spline.interpolate(x,y);
        derivativeFunc = dspsfunc.polynomialSplineDerivative();
        pc = getC();
    }

    /**
     * 计算
     * @param x
     * @return
     */
    public double value(double x){
        if (x<knots[0]){
            return pc[0][0]+pc[0][1]*(x-knots[0])+pc[0][2]*(x-knots[0])*(x-knots[0])+pc[0][3]*(x-knots[0])*(x-knots[0])*(x-knots[0]);
        }
        if (x>knots[knots.length-1]){
            return pc[knots.length-2][0]+pc[knots.length-2][1]*(x-knots[knots.length-2])+pc[knots.length-2][2]*(x-knots[knots.length-2])*(x-knots[knots.length-2])+pc[knots.length-2][3]*(x-knots[knots.length-2])*(x-knots[knots.length-2])*(x-knots[knots.length-2]);
        }
        return dspsfunc.value(x);
    }

    /**
     * 导数计算
     * @param x
     * @return
     */
    public double dvalue(double x){
        if (x<knots[0]){
            return pc[0][1]+2.0*pc[0][2]*(x-knots[0])+3.0*pc[0][3]*(x-knots[0])*(x-knots[0]);
        }
        if (x>knots[knots.length-1]){
            return pc[knots.length-2][1]+2.0*pc[knots.length-2][2]*(x-knots[knots.length-2])+3.0*pc[knots.length-2][3]*(x-knots[knots.length-2])*(x-knots[knots.length-2]);
        }
        return derivativeFunc.value(x);
    }


    /**
     * 返回各段的方程系数
     * @return
     */
    public double[][] getC(){
        PolynomialFunction[] ps = dspsfunc.getPolynomials();
        double[][] c = new double[ps.length][];
        for(int i=0;i<ps.length;i++){
            c[i] = ps[i].getCoefficients();
        }

        return c;
    }

    /**
     * 返回各段的端点值
     * @return
     */
    public double[] getKnots(){
        return dspsfunc.getKnots();
    }
}
